library(cjoint);library(haven);library(ggplot2);library(car)

##Read In Data
sppc_data <- read.qualtrics("YouGovCleaned.csv", 
                                responses=c("choice_1","choice_2","choice_3","choice_4","choice_5",
                                            "choice_6","choice_7","choice_8","choice_9","choice_10"), 
                                covariates=c("employ","joblooking","willingmove","birthyr","gender",
                                             "educ","faminc_new","pid7","ideo5","newsint"), respondentID="caseid",new.format=FALSE)

#Covariates need to be in the order they appear in the dataset. Commented is number of categories
table(sppc_data$pid7) #8
table(sppc_data$employ) #9
table(sppc_data$joblooking) #3
table(sppc_data$ideo5) #6
table(sppc_data$newsint) #1-4, 7
table(sppc_data$educ) #6
table(sppc_data$faminc_new) 

#Recode Data for Analysis
sppc_data$Education <- -1
sppc_data$Education[sppc_data$educ<4] <- 0
sppc_data$Education[sppc_data$educ>3] <- 1
table(sppc_data$Education,sppc_data$educ)
sppc_data$Education <- factor(sppc_data$Education,labels=c("No College Degree","College Degree"),levels=c(0,1))
prop.table(table(sppc_data$Education))

sppc_data$PID <- NA
sppc_data$PID[sppc_data$pid7==1] <- 1
sppc_data$PID[sppc_data$pid7==2] <- 1
sppc_data$PID[sppc_data$pid7==3] <- 1
sppc_data$PID[sppc_data$pid7==4] <- 2
sppc_data$PID[sppc_data$pid7==5] <- 3
sppc_data$PID[sppc_data$pid7==6] <- 3
sppc_data$PID[sppc_data$pid7==7] <- 3
table(sppc_data$PID)
sppc_data$PID <- factor(sppc_data$PID,labels=c("Democrat","Independent","Republican"),levels=c(1,2,3))
prop.table(table(sppc_data$PID))

sppc_data$Female <- -1
sppc_data$Female[sppc_data$gender==1] <- 0
sppc_data$Female[sppc_data$gender==2] <- 1
table(sppc_data$Female)
sppc_data$Female <- factor(sppc_data$Female,labels=c("Male","Female"),levels=c(0,1))
prop.table(table(sppc_data$Female))


sppc_data$Income <- -1
sppc_data$Income[as.numeric(sppc_data$faminc_new)<17] <- 4
sppc_data$Income[as.numeric(sppc_data$faminc_new)<10] <- 3
sppc_data$Income[as.numeric(sppc_data$faminc_new)<8] <- 2
sppc_data$Income[as.numeric(sppc_data$faminc_new)<5] <- 1
sppc_data$Income[as.numeric(sppc_data$faminc_new)==97] <- NA
table(sppc_data$Income)
sppc_data$Income <- factor(sppc_data$Income,labels=c("1st Quartile","2nd Quartile","3rd Quartile","4th Quartile"),levels=c(1,2,3,4))
table(sppc_data$Income)

sppc_data <- na.omit(sppc_data)

table(sppc_data$jantemp)
sppc_data$jantemp <- trimws(sppc_data$jantemp)
sppc_data$Average.January.Temperature <- factor(sppc_data$jantemp,levels=c("20 degrees Fahrenheit","30 degrees Fahrenheit","50 degrees Fahrenheit","65 degrees Fahrenheit"))
table(sppc_data$Average.January.Temperature)

sppc_data$salary <- trimws(sppc_data$salary)
sppc_data$Salary <- factor(sppc_data$salary,levels=c("$40,000","$70,000","$100,000","$130,000"))
table(sppc_data$Salary)

sppc_data$homeprice <- trimws(sppc_data$homeprice)
sppc_data$Typical.Home.Price <- factor(sppc_data$homeprice,levels=c("$170,000","$340,000","$410,000","$680,000"))
table(sppc_data$Typical.Home.Price)

sppc_data$companysize <- trimws(sppc_data$companysize)
sppc_data$Company.Size <- factor(sppc_data$companysize,levels=c("10 employees","500 employees","5,000 employees","Over 50,000 employees"))
table(sppc_data$Company.Size)

sppc_data$location <- trimws(sppc_data$location)
sppc_data$Location <- factor(sppc_data$location,levels=c("Major metropolitan area","Mid-size city","Small college town","Rural area"))
table(sppc_data$Location)

sppc_data$Recent.State.Social.Policy <- factor(sppc_data$statesocial, 
                              levels=c("Forming a commission to increase public understanding of the state's history and culture",
                                       "Affirming the ability of same-sex couples to adopt children",
                                       "Prohibiting same-sex couples from adopting children",
                                       "Banning the concealed carry of firearms",
                                       "Permitting the concealed carry of firearms",
                                       "Banning the teaching of certain ideas related to race and racial issues in the U.S. in high school history courses",
                                       "Requiring the teaching of race and racial issues in the U.S. in high school history courses",
                                       "Limiting abortion to the first 6 weeks of pregnancy with no exception for rape, incest, or health of the mother",
                                       "Reaffirming the right to abortion in the first trimester (and later if the health of the mother is in jeopardy)",
                                       "Allowing transgender student athletes to play on sports teams that match their gender identity, rather than the sex they were assigned at birth",
                                       "Requiring transgender student athletes to play on sports teams that match the sex they were assigned at birth, even if it does not match their gender identity"),
                              labels=c("Form Commission","Affirm LGB Adoption","Prohibit LGB Adoption","Ban Concealed Carry","Permit Concealed Carry",
                                       "Ban Teaching of Race","Require Teaching of Race","Limit Abortion","Affirm Right to Abortion","Let Trans Athletes Choose","Prohibit Trans Athletes Choose"))

sppc_data$Company.Culture <- factor(sppc_data$culture, 
                              levels=c("Employees are given many opportunities for advancement within the organization",
                              "The company seeks to provide employees with constructive feedback to foster their career growth",
                              "You will have many opportunities to collaborate with talented people",
                              "You will have the ability to work on a variety of tasks and develop your skills in many areas"),
                              labels=c("Advancement","Feedback","Collaborate","Task Variety"))

sppc_data$Recent.State.Economic.Policy <- factor(sppc_data$stateeconomic, 
                                                 levels=c("Forming a commission to study ways to strengthen economic growth",
                                                          "Decreasing K-12 per pupil spending from $12,000 to $6,000",
                                                          "Increasing K-12 per pupil spending from $12,000 to $18,000",
                                                          "Eliminating the state income tax",
                                                          "Increasing the income tax on millionaires",
                                                          "Expanding Medicaid coverage to 140% of the poverty level",
                                                          "Keeping Medicaid coverage at 100% of the poverty level rather than expanding it to 140% of the poverty level",
                                                          "Keeping the minimum wage at $7.25 per hour rather than raising it to $15 per hour", 
                                                          "Raising the minimum wage from $7.25 to $15 per hour"),
                                                 labels=c("Form Commission","Decrease K-12 Spending","Increase K-12 Spending","Eliminate Income Tax","Tax Millionaires",
                                                          "Expand Medicaid","Stable Medicaid",
                                                          "Keep Minimum Wage","Raise Minimum Wage"))
table(sppc_data$Recent.State.Economic.Policy)

sppc_data$Presidential.Election.Results <- factor(sppc_data$election, 
                                    levels=c("In a state that voted heavily for Donald Trump",
                                             "In a state that Donald Trump barely won",
                                             "In a state that Joe Biden barely won",
                                             "In a state that voted heavily for Joe Biden"),
                                    labels=c("Strong Trump","Weak Trump","Weak Biden","Strong Biden"))
table(sppc_data$Presidential.Election.Results)

conjoint_theme <- theme(panel.border = element_rect(colour = "black",fill=NA,size=1.25),
                        panel.background = element_rect(fill = "white",
                                                        colour = "white",
                                                        size = 0.5, linetype = "solid"),
                        panel.grid.major = element_line(size = 0.5, linetype = 'solid',
                                                        colour = "grey80"), 
                        panel.grid.minor = element_line(size = 0.25, linetype = 'solid',
                                                        colour = "grey80"),
                        legend.position="none",
                        strip.background = element_rect(colour="black", fill=NA,size=1.25),
                        strip.text.x = element_text(
                          size = 14, color = "black", face = "bold"),
                        axis.text.y=element_text(hjust=0),
                        axis.text=element_text(colour="black"),
                        axis.title = element_text(size = 14))

sppc_data$RecentStateSocialPolicy <- sppc_data$Recent.State.Social.Policy
sppc_data$RecentStateEconomicPolicy <- sppc_data$Recent.State.Economic.Policy
sppc_data$TypicalHomePrice <- sppc_data$Typical.Home.Price
sppc_data$AverageJanuaryTemperature <- sppc_data$Average.January.Temperature
sppc_data$CompanyCulture <- sppc_data$Company.Culture
sppc_data$CompanySize <- sppc_data$Company.Size
sppc_data$PresidentialElectionResults <- sppc_data$Presidential.Election.Results


##Figure 1
results <- amce(selected ~ RecentStateSocialPolicy+RecentStateEconomicPolicy+Salary+TypicalHomePrice+AverageJanuaryTemperature+CompanyCulture+CompanySize+Location+PresidentialElectionResults, 
                data=sppc_data,respondent.id="respondent")

summary(results)
pdf("Figures/YG_Selected.pdf",width=8,height=11,paper='special')
plot(results,plot.theme=conjoint_theme,group.order=c("RecentStateSocialPolicy",
                                                     "RecentStateEconomicPolicy",
                                                     "Salary",
                                                     "TypicalHomePrice",
                                                     "AverageJanuaryTemperature",
                                                     "CompanyCulture",
                                                     "CompanySize",
                                                     "Location",
                                                     "PresidentialElectionResults"),
     xlab="Change in Pr(Selected)")
dev.off()


##Figure 2
round(coeftest(lm(selected ~ Female*Average.January.Temperature+Female*Company.Culture+Female*Company.Size+Female*Location+Female*Presidential.Election.Results+Female*Recent.State.Economic.Policy+Female*Recent.State.Social.Policy+Female*Salary+Female*Typical.Home.Price, 
                  data=sppc_data),
               vcov = vcovCL,
               type = "HC1",
               cluster = ~respondent),2)


pdf("SocialConjointYGResultsGender.pdf",width=6,height=6,paper='special') 
par(mar=c(3,10,2,1),mfrow=c(1,1))
plot(NA,xlim=c(-.3,.2),ylim=c(1,5),yaxt="n",ylab="",xlab="Difference in in Pr(Selected)\n",
     main="Gender")
points(x=c(-0.01,0.07,-0.01,-0.09,-0.06),
       y=c(5:1),pch=19)
segments(x0=-0.01+1.96*0.04,
         x1=-0.01-1.96*0.04,y0=5,y1=5,lwd=2)
segments(x0=0.07+1.96*0.04,
         x1=0.07-1.96*0.04,y0=4,y1=4,lwd=2)
segments(x0=-0.01+1.96*0.04,
         x1=-0.01-1.96*0.04,y0=3,y1=3,lwd=2)
segments(x0=-0.09+1.96*0.04,
         x1=-0.09-1.96*0.04,y0=2,y1=2,lwd=2)
segments(x0=-0.06+1.96*0.04,
         x1=-0.06-1.96*0.04,y0=1,y1=1,lwd=2)
abline(v=0,col="lightgray",lty=2,lwd=3)
axis(2,5:1,labels=c("Prohibit LGB Adoption","Ban Concealed Carry","Ban Teaching Race",
                    "Limit Abortion","Prohibit Trans Athletes"),las=1)
dev.off()

round(coeftest(lm(selected ~ Education*Average.January.Temperature+Education*Company.Culture+Education*Company.Size+Education*Location+Education*Presidential.Election.Results+Education*Recent.State.Economic.Policy+Education*Recent.State.Social.Policy+Education*Salary+Education*Typical.Home.Price, 
                  data=sppc_data),
               vcov = vcovCL,
               type = "HC1",
               cluster = ~respondent),2)
pdf("SocialConjointYGResultsEduc.pdf",width=6,height=6,paper='special') 
par(mar=c(3,10,2,1),mfrow=c(1,1))
plot(NA,xlim=c(-.3,.2),ylim=c(1,5),yaxt="n",ylab="",xlab="Difference in Pr(Selected)\n",
     main="College Degree")
points(x=c(0.01,0.06,0.05,0.03,-0.02),
       y=c(5:1),pch=19)
segments(x0=0.01+1.96*0.04,
         x1=0.01-1.96*0.04,y0=5,y1=5,lwd=2)
segments(x0=0.06+1.96*0.04,
         x1=0.06-1.96*0.04,y0=4,y1=4,lwd=2)
segments(x0=0.05+1.96*0.04,
         x1=0.05-1.96*0.04,y0=3,y1=3,lwd=2)
segments(x0=0.03+1.96*0.04,
         x1=0.03-1.96*0.04,y0=2,y1=2,lwd=2)
segments(x0=-0.02+1.96*0.04,
         x1=-0.02-1.96*0.04,y0=1,y1=1,lwd=2)
abline(v=0,col="lightgray",lty=2,lwd=3)
axis(2,5:1,labels=c("Prohibit LGB Adoption","Ban Concealed Carry","Ban Teaching Race",
                    "Limit Abortion","Prohibit Trans Athletes"),las=1)
dev.off()

round(coeftest(lm(selected ~ Income*Average.January.Temperature+Income*Company.Culture+Income*Company.Size+Income*Location+Income*Presidential.Election.Results+Income*Recent.State.Economic.Policy+Income*Recent.State.Social.Policy+Income*Salary+Income*Typical.Home.Price, 
                  data=sppc_data),
               vcov = vcovCL,
               type = "HC1",
               cluster = ~respondent),2)
pdf("SocialConjointYGResultsIncome.pdf",width=6,height=6,paper='special') 
par(mar=c(3,10,2,1),mfrow=c(1,1))
plot(NA,xlim=c(-.3,.2),ylim=c(1,5),yaxt="n",ylab="",xlab="Difference in Pr(Selected)\n",
     main="High Income")
points(x=c(-0.06,0.04,0,0.02,-0.02),
       y=c(5:1),pch=19)
segments(x0=-0.06+1.96*0.05,
         x1=-0.06-1.96*0.05,y0=5,y1=5,lwd=2)
segments(x0=0.04+1.96*0.06,
         x1=0.04-1.96*0.06,y0=4,y1=4,lwd=2)
segments(x0=0+1.96*0.05,
         x1=0-1.96*0.05,y0=3,y1=3,lwd=2)
segments(x0=0.02+1.96*0.05,
         x1=0.02-1.96*0.05,y0=2,y1=2,lwd=2)
segments(x0=-0.02+1.96*0.02,
         x1=-0.02-1.96*0.02,y0=1,y1=1,lwd=2)
abline(v=0,col="lightgray",lty=2,lwd=3)
axis(2,5:1,labels=c("Prohibit LGB Adoption","Ban Concealed Carry","Ban Teaching Race",
                    "Limit Abortion","Prohibit Trans Athletes"),las=1)
dev.off()








